## Replication for Unilateral Powers, Public Opinion, and the Presidency by
## Andrew Reeves and Jon Rogowski Use this code to replication column 1 of
## Tables 1 and 2, Table 3 in the main text, and column 1 of Tables 1, 2,
## and 3 of the Online SI.

rm(list = ls())
library(psych)
library(apsrtable)

## We load the full CCES data file and subset out the variables used for this analysis
#ccesData <- read.spss("data/cces/WAS/CCES13_WAS/CCES13_WAS_OUTPUT.sav", to.data.frame=TRUE)
#ccesData <- ccesData[,c("caseid", "CC312a", "CC334A", "pid7", "race", "birthyr", "gender", "educ",
#                        "faminc", "WAS023", "WAS024", "WAS025", "WAS026", "WAS009", "WAS008",
#                        "WAS010", "weight")]
#save(ccesData, file = "rr_cces.Rdata")
load("rr_cces.Rdata")
## Recode variables for analyses
## Obama approval
ccesData$approval <- as.numeric(ccesData$CC312a)
ccesData$approval <-
  recode(ccesData$approval, "5=NA;1=4;2=3;3=2;4=1")
## Ideology
ccesData$ideology <-
  as.numeric(
    recode(
      ccesData$CC334A, "'Very Liberal'=7;'Liberal'=6;'Somewhat Liberal'=5; 'Middle of the Road'=4;'Somewhat Conservative'=3; 'Conservative'=2;'Very Conservative'=1;'Not sure'=NA"
    )
  )
## Party id
ccesData$dem <- rep(0, nrow(ccesData))
ccesData$dem[ccesData$pid7 %in% c("Strong Democrat", "Not very strong Democrat", "Lean Democrat")] <-
  1
ccesData$gop <- rep(0, nrow(ccesData))
ccesData$gop[ccesData$pid7 %in% c("Strong Republican", "Not very strong Republican", "Lean Republican")] <-
  1
ccesData$pid <-
  as.numeric(
    recode(
      ccesData$pid7,"'Strong Republican'=1;'Not very strong Republican'=2;'Lean Republican'=3;'Strong Democrat'=7;'Not very strong Democrat'=6; 'Lean Democrat'=5;'Independent'=4;'Not sure'=NA"
    )
  )
## Race
ccesData$black <- rep(0, nrow(ccesData))
ccesData$black[ccesData$race == "Black"] <- 1
ccesData$hispanic <- rep(0, nrow(ccesData))
ccesData$hispanic[ccesData$race == "Hispanic"] <- 1
## Age
ccesData$age <- 2013 - ccesData$birthyr
ccesData$senior <- rep(0, nrow(ccesData))
ccesData$senior[ccesData$age > 64] <- 1
## Gender
ccesData$male <- rep(0, nrow(ccesData))
ccesData$male[ccesData$gender == "Male"] <- 1
## BA Education
ccesData$education <- as.numeric(ccesData$educ)
## Income
ccesData$income <- as.numeric(ccesData$faminc)
ccesData$income[ccesData$income == 19] <- NA
## Rule of Law
## A function to recode variables
flipScale <- function(x) {
  tmp <- as.numeric(x)
  tmp <- recode(tmp, "1=4;2=3;3=2;4=1")
  return(tmp)
}
ccesData$law.scale <-
  alpha(data.frame(
    flipScale(ccesData$WAS023), flipScale(ccesData$WAS024), flipScale(ccesData$WAS025), flipScale(ccesData$WAS026)
  ))$score

## Recode our dependant variables.  Text of questions included.  Variables with .bin suffixes are binary indicators of support.
## The president should be able to appoint judges of his choosing regardless of
## whether the U.S. Senate agrees with his selections.
ccesData$judges.bin <- rep(0, nrow(ccesData))
ccesData$judges.bin[as.numeric(ccesData$WAS008) < 3] <- 1
ccesData$judges <-
  as.numeric(
    recode(
      ccesData$WAS008, "'Strongly Agree'=4;'Somewhat Agree'=3;'Somewhat Disagree'=2;'Strongly Disagree'=1"
    )
  )

## A president should have the authority to decide
## how executive agencies will implement bills passed by Congress.
ccesData$agencies.bin <- rep(0, nrow(ccesData))
ccesData$agencies.bin[as.numeric(ccesData$WAS009) < 3] <- 1
ccesData$agencies <-
  as.numeric(
    recode(
      ccesData$WAS009, "'Strongly Agree'=4;'Somewhat Agree'=3;'Somewhat Disagree'=2;'Strongly Disagree'=1"
    )
  )

## A president should have the right to enact policies without having those
## policies voted on by Congress.
ccesData$execaction.bin <- rep(0, nrow(ccesData))
ccesData$execaction.bin[as.numeric(ccesData$WAS010) < 3] <- 1
ccesData$execaction <-
  as.numeric(
    recode(
      ccesData$WAS010, "'Strongly Agree'=4;'Somewhat Agree'=3;'Somewhat Disagree'=2;'Strongly Disagree'=1"
    )
  )
## Replicate results from column 1 of Table 1
weighted.mean(ccesData$execaction.bin, ccesData$weight)
weighted.mean(ccesData$agencies.bin, ccesData$weight)
weighted.mean(ccesData$judges.bin, ccesData$weight)
weighted.mean(recode(ccesData$approval, "1=0;2=0;3=1;4=1"), ccesData$weight, na.rm = TRUE)
## Replicate results from column 1 of Table 2
weighted.mean(ccesData$execaction.bin[ccesData$approval > 2], ccesData$weight[ccesData$approval > 2], na.rm = TRUE)
weighted.mean(ccesData$execaction.bin[ccesData$approval < 3], ccesData$weight[ccesData$approval < 3], na.rm = TRUE)
weighted.mean(ccesData$judges.bin[ccesData$approval > 2], ccesData$weight[ccesData$approval > 2], na.rm = TRUE)
weighted.mean(ccesData$judges.bin[ccesData$approval < 3], ccesData$weight[ccesData$approval < 3], na.rm = TRUE)
weighted.mean(ccesData$agencies.bin[ccesData$approval > 2], ccesData$weight[ccesData$approval > 2], na.rm = TRUE)
weighted.mean(ccesData$agencies.bin[ccesData$approval < 3], ccesData$weight[ccesData$approval < 3], na.rm = TRUE)


## Replicate Table 3 (note that these results are also column 1 of Table 2.1, 2.2, and 2.3 in the Online SI)
## Specify the model
themodel.cong <- "execaction.bin ~ law.scale + approval  + ideology"
themodel.judges <- "judges.bin ~ law.scale + approval + ideology"
themodel.agencies <-
  "agencies.bin ~ law.scale + approval + ideology"
## Estimate the model
cces.model.cong <- glm(
  as.formula(themodel.cong), data = ccesData,
  family = binomial(link = "logit"),
  weights = weight
)
cces.model.judges <-
  glm(
    as.formula(themodel.judges), data = ccesData,
    family = binomial(link = "logit"),
    weights = weight
  )
cces.model.agencies <-
  glm(
    as.formula(themodel.agencies), data = ccesData,
    family = binomial(link = "logit"),
    weights = weight
  )
## Make Table 3
table3 <-
  apsrtable(
    cces.model.cong, cces.model.agencies, cces.model.judges,
    model.names = c(
      "Unilateral Policymaking","Bureaucratic Implementation", "Judicial Appointments"
    ),
    coef.names = c(
      "Intercept", "Belief in rule of law", "Presidential Approval",
      "Ideology"
    ),
    caption.position = "below",
    caption = "\\textbf{A Model of Support for Unilateral Actions: Unilateral Policymaking,
    Bureaucratic Implementation, and Judicial Appointments.} Each columns show results for support
    of a different type of unilateral power.  Both presidential approval and beliefs in rule of law
    are significantly related to attitudes toward executive action.  The dependent variable is whether
    a respondent agrees that, ``the president should be able to appoint judges of his choosing regardless
    of whether the U.S. Senate agrees with his selections.''",
    label = "t:execaction",
    Sweave = TRUE
  )

table3

## Make Figure 1 in main text
## be sure to download this hack from the dataverse
source("plot.ci2.R")
detach("package:psych", unload=TRUE)
makePIDfigures <- function(dv, ivs = c("approval", "ideology", "law.scale"), data = ccesData, ...){
  data <- data[,c(dv, ivs, "weight")]
  data <- subset(data, complete.cases(data))
  z.out <- zelig(as.formula(paste(dv, "~", paste(ivs, collapse = "+"))),
                 model = "logit",
                 data = data,
                 weights = data$weight)
  x.out <- setx(z.out, approval = seq(from = 1, to = 4, by = .25))
  s.out <- sim(z.out, x = x.out)
  plot.ci2(s.out, 
           ci=c(95), 
           ylim = c(0,1),
           legpos = c(200,200),
           ...)
  par(las=1)
  axis(side = 1, at = c(1,2,3,4), labels = c("","","",""))
  mtext("Strongly\n Approve", 1, at = 4, line = 1, cex = .5)
  mtext("Strongly\n Dispprove", 1, at = 1, line = 1, cex = .5)
}
par(mfrow = c(2,3), 
    oma = c(2, 8.25, 2, 2),
    mar = c(6, 2, 0, 0), 
    mgp = c(2, 1, 0), 
    xpd = NA)
makePIDfigures(dv = "execaction.bin", main = "Unilateral Policymaking", ylab = "Probability of Support", xlab = "Approval")
mtext("Obama\nApproval", 2, line = 4, cex = .75)
par(las=1)
makePIDfigures(dv = "judges.bin", thesurvey.month = "Nov13", ylab = "", xlab = "Approval", main = "Judicial Appointments")
makePIDfigures(dv = "agencies.bin", thesurvey.month = "Nov13", ylab = "", xlab = "Approval", main = "Bureaucratic Implementation")
## Rule of Law
makeROLfigures <- function(dv, ivs = c("approval", "ideology", "law.scale"), data = ccesData, ...){
  data <- data[,c(dv, ivs, "weight")]
  data <- subset(data, complete.cases(data))
  z.out <- zelig(as.formula(paste(dv, "~", paste(ivs, collapse = "+"))),
                 model = "logit",
                 data = data,
                 weights = data$weight)
  x.out <- setx(z.out, law.scale = seq(from = 1, to = 4, by = .25))#, approval = 2.5)
  s.out <- sim(z.out, x = x.out)
  plot.ci2(s.out, 
           ci=c(95), 
           ylim = c(0,1),
           legpos = c(200,200),
           ...)
  axis(side = 1, at = c(1,2,3,4), labels = c("1","2","3","4"))  
}

makeROLfigures(dv = "execaction.bin", thesurvey.month = "Nov13", ylab = "Probability of Support", xlab = "Belief in the Rule of Law")
mtext("Belief in\nRule of Law", 2, line = 4, cex = .75)
par(las=1)
makeROLfigures(dv = "judges.bin", thesurvey.month = "Nov13", ylab = "", xlab = "Belief in the Rule of Law")
makeROLfigures(dv = "agencies.bin", thesurvey.month = "Nov13", ylab = "", xlab = "Belief in the Rule of Law")
dev.off()
